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The numerical extraction of resonant states of open quantum systems is usually a difficult problem. 
Regularization techniques, such as the mapping to complex coordinates or the addition of Complex 
Absorbing Potentials are typically employed, as they render resonant wavefunctions localized and 
therefore normalizable. Physically relevant metastable states have energies that do not depend on 
the chosen regularization method. Their identification therefore involves cumbersome comparisons 
between multiple regularised calculations, often performed graphically, which require fine-tuning and 
specific intuition to avoid approximated, if not wrong, results. In this Letter, we define an operator 
that explicitly measures such invariance, valid for any arbitrary mapping of spatial coordinates. 
Resonant states of the system can eventually simply be identified evaluating the expectation value 
of this operator. Our method eases the extraction of resonant states even for numerical potentials 
that are difficult to scale to complex coordinates, and avoids the need for ad hoc complex absorbing 
potentials. We provide explicit evidence of our findings discussing one-dimensional case-studies, also 
in the presence of external electric fields. 


Resonant states are ubiquitous in Quantum Physics. 
Also referred to as “Gamow vectors” or “Siegert states”, 
they can be defined as solutions of the time-independent 
Schrodinger equation subject to outgoing boundary con¬ 
ditions. Described at first by Gamow [1] via quasi¬ 
stationary states, the concept of resonant states has 
been widely developed in the field of atomic and nuclear 
physics (see e.g. Ref. [2]), then adopted for the analysis 
of scattering properties of quantum systems with open 
boundaries [3]. Various literature has shown that Siegert 
states encode in compact form the response properties 
of a system [4, 5]. In particular, the analytic structure 
of the resolvent operator (i.e. the Green’s function) is 
completely determined by resonant energies and wave- 
functions. In the words of Ref. [6], resonant states ex¬ 
pansions offer the “possibility of a unified description of 
bound states, resonances, and continuum spectrum in 
terms of a purely discrete set of states”. For one-body 
Hamiltonians of quantum systems, the identification of 
their resonant energies and wavepackets is therefore of 
paramount importance. 

Resonant wavefunctions exhibit complex wavenumbers 
k = ha — iki and divergent asymptotes oc exp(ifc • x). 
As such, they are not to be found in the Hilbert space 
of square-integrable functions. Yet, they find rigor¬ 
ous mathematical foundations in the domain of non- 
Hermitian Quantum Mechanics [7]. In order to ease 
their numerical treatment, methods have been proposed, 
which render resonant states square-integrable, allowing 
their computation under bound-state-like (i.e. Dirichlet) 
boundary conditions. These regularization methods im¬ 
ply a modification of the original Hamiltonian of the sys¬ 
tem, driven by a set of continuous parameters (usually 
called r] or 9, see e.g. [8]), which eventually leads to a 
complex-valued operator, i.e. explicitly non-Hermitian. 
In this way, localized - thus square-integrable - eigen¬ 
vectors with complex eigenvalues may show up in the 


spectrum of the resulting operator. In this context, phys¬ 
ically relevant resonant states are those having an energy 
that is invariant with respect to the chosen regularization 
method. 

The identification of resonant energies requires com¬ 
parisons between multiple regularised calculations, typ¬ 
ically performed graphically (via so-called rj- or 6- 
trajectories). Such a comparison is generally difficult to 
be performed and requires highly precise calculations to 
identify stable points in the spectra of several complex 
Hamiltonians. This is especially true when the potential 
is only known numerically. Moreover, even when a stable 
point is found, little is known about the numerical qual¬ 
ity of the corresponding eigenvector, that may depend on 
computational parameters such as the size of simulation 
box or the choice of the numerical basis set. 

In this Letter, we revisit the properties of Siegert states 
under arbitrary parametric transformations of spatial co¬ 
ordinates. We eventually introduce an operator whose 
quantum expectation value is explicitly associated to the 
variation of the energy with respect to the parameters of 
the coordinate mapping. Such operator can be related 
to a generalisation of the classical Virial Theorem for 
stationary states. This procedure provides a reliable and 
rigorous approach to identify resonant states without the 
need neither of explicit variations of the parameters nor 
the analytic continuation of the numerical potential in 
the complex plane. The method we propose is crucial in 
indicating which states have some physical interest and, 
at the same time, provides an estimate of the accuracy of 
the computational treatment. By presenting some illus¬ 
trative examples involving one-dimensional models, we 
demonstrate that our approach can be used to single out 
the states that are more relevant in determining the lin¬ 
ear response of open quantum systems. 

In order to motivate the interest of our results, let 
us first illustrate the main advantages and drawbacks of 
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the most popular regularization techniques. In a num¬ 
ber of numerical investigations, resonant states are com¬ 
puted through the Complex Absorbing Potential (CAP) 
method [9], where a complex potential is added to the 
Hamiltonian such as to absorb the decaying particle de¬ 
scribed by the outgoing resonant state. Within this ap¬ 
proach, the original one-body potential is not modified, 
which makes this approach suitable for numerical poten¬ 
tials [10]. A resonant state is then identified as a CAP- 
independent state, and its energy is often found by ver¬ 
ifying numerically its invariance with respect to varia¬ 
tions of the CAP strength (the ry-trajectories). However, 
there is no unambiguous recipe for the CAP to ensure 
that some resonant states appear in the spectrum of the 
non-Hermitian Hamiltonians. The absorbing boundary 
described by the CAP may induce artificial “reflections” 
of the resonant state wavefunction at the boundaries of 
the simulation domain, thereby altering their energy as 
well as their shape within the quantum device. 

A method which is based on rigorous mathematical 
foundations is the well-known Complex Scaling Method 
(CSM) [11], in which spatial coordinates are “scaled” by 
a complex factor, x —> xe'®. All the resonant states for 
which ki/kR < tan0 become localized, and their ener¬ 
gies are ^-independent. Instead, continuum states show 
up along straight lines, rotated by an angle —29 with 
respect to the real axis. If the Hamiltonian has only 
one threshold energy {i.e. P(oo), see Ref. [12]), resonant 
energies can be identified by simply looking at their po¬ 
sition with respect to the rotated continuum. However, 
even for the CSM, the identification of a resonant energy 
often relies on its numerical invariance with respect to 6 
(0-trajectories). 

Despite its conceptual simplicity, CSM is unfit for the 
treatment of generic numerical potentials with large spa¬ 
tial extension, as it introduces high-frequency oscillations 
of the potential far from the fixed point of the scaling 
transformation. This is true even for analytic potentials. 
To illustrate this point, it is enough to consider the sim¬ 
plest prototype of a localized, smooth function, a Gaus¬ 
sian G{a,Xc',x), undergoing a complex scaling transfor¬ 
mation centered at the origin: 


G(cr, Xc', xe'^) = exp 


(xe'® — XcY 
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( 1 ) 


Such function can model a “diffusive center” placed at the 
position Xc in the simulation domain. When Xc/u ^ 1, 
as in spatially extended systems like electronic poten¬ 
tials with several diffusive centers, the complex scaling 
transformation induces high-frequency (as well as high- 
amplitude) oscillations on the potential, very difficult to 
be captured numerically. The potential becomes so oscil¬ 
lating that an accurate numerical treatment is unfeasible 
even for computational domains of moderate size: the 
numerical basis set should be able to capture the large 
and rapid oscillations both of the potential and of the 


eigenvectors, making the computational cost overwhelm¬ 
ing. 

An elegant generalization of the CSM exists, still based 
on rigorous foundations. Referred to as “Reflection-Free 
CAP” in Ref. [13], or as “Smooth Exterior complex Scal¬ 
ing Method” (SESM), this approach somehow couples 
the CAP method with the CSM. The SESM stems from 
the coordinate transformation 


X -> F^{x) , 


( 2 ) 


^ = 0 being the identity transformation. A rigorous, 
non-Hermitian Hamiltonian can be obtained out of the 
reparametrization F^(x), from which resonant states can 
be extracted. The function F^{x) is generally chosen to 
tend asymptotically to x e‘^ when x —> ±oo, thereby rec¬ 
onciling with the CSM. When the family of functions 
can be chosen so that F^{x) = x where V(x) yf 0, the 
potential of the SESM Hamiltonian can be left unsealed. 
Also non-local {e.g. many-body) potentials, whose ana¬ 
lytic continuation to complex coordinates might be cum¬ 
bersome, can be treated with this method. As in the 
CSM, resonant state energies have to be ^-independent. 
However, as in the CAP method, continuum states can¬ 
not be easlily excluded, hence resonant energies have to 
be found by explicitly verifying their independence with 
respect to the ^ parameter space (typically, only 9 is con¬ 
sidered [14]). 

For multi-centered potentials, methods like CAP or 
SESM seem very interesting, as spatial coordinates can 
be left unsealed in the inner region of the simulation do¬ 
main, whose extension is related to that of the potential. 
In particular, the bound states whose support is con¬ 
tained in the unsealed region are not modified. This fact 
has a remarkable practical consequence: bound states 
of Hamiltonian can be first extracted with ^ = 0, then 
regarded as exact eigenstates of the SESM Hamiltonian, 
provided that their support is within the unsealed region. 

Aside from the SESM or CSM, let us now consider the 
coordinate mapping of Eq. (2) on a completely general 
ground, by assuming a generic form of the function F^. 
In what follows we adopt the notation 




= (3) 


For the normalization of bound states to be pre¬ 
served, wavefunctions have to transform as tp ^ ip^ = 
^yy^ip{F^). Such transformation induces the following 
modification on the position and momentum operators, 
X and P = —idx = ife (cf. Ref. [13]): 


A 


P 


= A[P^] = F^{X) , 

A = HM = 


( 4 ) 


K- 


vrrvk' 


Since the transformation (2) has to preserve phys¬ 
ical quantities, physically meaningful states are ex- 
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pected to be eigenstates of the transformed Hamiltonian a dilation of the original Hamiltonian F\{x) = e^x 
= T[/j] + V(F^), such that We have 


(V’dV'c) (dddc) 


(5) 


dH. 


dX 


-T-ltt + xV' = 


x,P =i r 


W,H 


( 10 ) 


A=0 


where T = 12, and are the right and left 

eigenvectors respectively, and the last equality derives 
from the Hellmann-Feynman theorem. Eq. (5) is of 
course also valid for physical eigenstates with complex 
energy. By expressing 


tj 

(dd^ldc) = 


( 6 ) 


du 


/?(w) (dd 


sJM 




we can identify an operator whose expectation value has 
to be zero on states having an energy that is invariant 
under reparametrizations (2): 


^ = ;7[F4,/d + i/i[/d^+^C>ii[/d^+t>2[/d5^ (7) 


where 


[/[F,/] = E'(F)F+i/|^-i^| , (8) 

Iff If If 

Uilf] ^ ^ "2:^ ’ ^ 2^ ■ 


All terms can be evaluated in closed-form, except for 
y'(E’^). However, in practical applications of the SESM, 
it should not be necessary to evaluate V', as F^ is de¬ 
signed such that yf 0 only where V = 0. 

The scaling behaviour of wavefunctions suggests that 
it exists another operator enabling us to single out phys¬ 
ical states. Since stems from the transformation of 
X and P in Eq. (4), a physically meaningful state 
has to be compatible with the canonical commutation 
relation of Quantum Mechanics. In other terms, for the 
(regularized) normalization (dcIV'?) to be consistent, the 


relation 




= i (V'jl'dc) should hold within 


the chosen regularization scheme. Thus, a physical eigen¬ 
state must also satisfy 


(dd^ddc) 

(ddd?) 


Oi = l-FA 
h 


+ (9) 

k 


where W = \ (^P + PX^ = — ^ (^x^ — ^x^ is the 
Weyl-quantized form of the dilation generator x • p. The 
classical virial theorem shows that the latter quantity 
is conserved on stationary orbits. The second member, 
which has to be zero on physical states, corresponds to 
the operator already presented in Ref. [12], called for this 
reason Complex Virial Operator. In a numerical compu¬ 
tation, its expectation value is related to the “pressure” 
exerted by the state on the boundaries of the simulation 
domain. Being a commutator with H, the last mem¬ 
ber of Eq. (10) is of course zero on any eigenstate of H. 
However, Eq. (10) only holds when evaluated on states 
satisfying Eq. (9): on such states the operator TV is a 
truly conserved quantity. This condition can be easily 
generalized to arbitrary transformations: for any choice 
of F^, each eigenstate of satisfying Eq. (9) will have 
a ^-independent energy. 

The operators in Eqs. (7, 9) constitute the main re¬ 
sults of this paper. The quantity {'tp^\0^\'ip^) is indirectly 
associated with the variation of the energy {'tjj^\F[^\'tp^) in¬ 
duced by the coordinate mapping (2). Eq. (9) can there¬ 
fore be used as an alternative to Eq. (5), with the advan¬ 
tage that no explicit derivative with respect to any of the 
^ parameters is needed. In both cases, the fulfillment of 
the equation provides an actual criterion for distinguish¬ 
ing, a posteriori, Siegert states from continuum states. 

As physical results should not depend on computa¬ 
tional parameters such as the simulation box size, the 
degree of fullfillment of Eqs. (5, 9) in numerical treat¬ 
ments provides an indication of the quality of the compu¬ 
tational setup. Matricial representations of these equa¬ 
tions can be easily evaluated in any basis set; we point 
out, however, that compositions of operators have to be 
done explicitly before projection onto the basis set. 

In order to demonstrate the significance of the infor¬ 
mation obtained from Eq. (7), we first discuss its numer¬ 
ical application to one-dimensional potentials, made up 
by superpositions of Gaussian functions as indicated in 
(1). Being v{xc,x) = [G(4, -Xc,x) + G(4,Xc;x)], we 

define (see Fig. 1-a) 

Vs{x) = v{45, x) , Vc{x) = v{20,x) — 0.2Vs{x). (11) 


where 


1-fi 




If has no boundary terms. 


the above condition is evidently satisfied and no explicit 
regularization is needed. It should be noted that, even for 
^ = 0, no regularization is possible for continuum states 
to satisfy Eq. (9). 

Moreover, it can be shown that, for any eigenstate of 
FI^, Eq. (9) implies Eq. (5). For instance, let us consider 


The Hamiltonian with 14 models two diffusive centers far 
from each other, yet close enough to exhibit a non-trivial 
scattering structure. The other toy potential 14, tailored 
to preserve the degeneracy of the bound eigenspaces, 
models a charge transfer from a “core” region, to the 
“shell” potential 14. These two systems therefore present 
very similar bound states but rather different spectra of 
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resonant states. The one-dimensional Hamiltonians are 
discretized in a Daubechies wavelet basis set as described 
in Ref. [12]. The SESM method is applied [15]. 

Fig. 1-b provides the evidence that bound and resonant 
states are stable with respect to variations of the scaling 
parameter 9. This fact is clearly confirmed by the corre¬ 
sponding value of obtained from Eq. (7). In Fig. 1-c, 
we show that for the potential Vs the SESM outperforms 
the CSM. Since the Gaussians are far from each other, 
in case of the CSM the potential and the eigenfunctions 
oscillate so strongly that very small grid spacings are 
required for the reliable identification of resonant ener¬ 
gies. Even with grids five times denser than the one used 
at 0 = 0 - already dense enough to correctly discretize 
Vs{xe'^) - the CSM still provides incorrect estimates of 
resonant energies. This fact is further confirmed by the 
corresponding values of the complex virials. 

As mentioned in the introduction, resonant states can 
be used to obtain a compact representation of the Green’s 
function. In systems with open boundaries, such a com¬ 
pact representation is very useful to express linear re¬ 
sponse functions and the derived optical spectroscopic 
properties in an optimal way. A good indicator of the 
quality of a discrete basis set for optical properties is 
the fulfillment of the Thomas-Reiche-Kuhn sum rule (or 
/-sum rule, see e.g. Ref. [17]), which relates the first mo¬ 
mentum of the oscillator strengths to that of the equilib¬ 
rium density of the system, i.e. the number of states be¬ 
low the Fermi level. The information provided by Eq. (7) 
is of great utility in this case: we have plotted in Fig. 1 the 
fulfillment of the /-sum rule as a function of the number 
of unoccupied states considered, ordered by for the 
Hamiltonian with the 14 potential. With such ordering, 
we are able to identify the minimal set of states allowing 
the fulfillment of the sum rule up to an excellent accuracy, 
independently of the simulation box size. This cannot be 
achieved with the low-energy eigenstates of the original, 
unsealed Hamiltonian, which suffer from the well-known 
continuum collapse [18]. Equivalent results can also be 
obtained by using the operator of Eq. (9). 

Our approach can also be used to identify field-induced 
metastable states of Hamiltonians in non-trivial environ¬ 
ments. A notable example of such metastability is given 
by the Stark ionization of molecules, recently studied 
within the framework of Density Functional Resonance 
Theory [19]. The reliable determination of metastable 
states that have to be occupied is admittedly a problem 
in the latter treatment. However, our method provides 
a natural solution to this problem. As an example, in 
Fig. 2, we show the values of (Oj) of Eq. (9) for a model 
Hamiltonian with a “soft-Coulomb” ID potential under 
external electric field of intensity 8, see e.g. Ref. [16]. 
States with smallest [(O^)] are those originating from a 
perturbation of the Rydberg states at £" = 0, therefore 
not belonging to the set of “continuum” states of the 
SESM Hamiltonian. This approach would be very help- 
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FIG. 1. (a) ID model potentials of Eq. (11). (b) Spectra 
obtained with different SESM tranformation, compared with 
the corresponding value of ^ from Eq. (7). The simulation 
box considered has size L = 300 AU. (c) SESM is compared 
against CSM for the Hamiltonian with 14 potential. The 
error in the identification of some resonant energies (pointed 
by the arrows in panel (b)) is plotted against the number N 
of degrees of freedom used in the simulation. Reference values 
are extracted with the SESM at A = 3000. (d) Fulfillment 
of the /-sum rule for the 14 potential as a function of the 
number of empty states considered for different L, ordered by 
I ^ I. Fermi level is chosen such as to have 8 occupied states. 

ful in identifying the physical metastable channels of open 
quantum devices under external electric fields [10]. 

To summarize, we have presented a simple and ro¬ 
bust method to numerically identify the Siegert states 
of open quantum systems undergoing a generic coordi¬ 
nate reparametrization. Such method, based on rigorous 
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FIG. 2. Measurement of (O 5 ) from Eq. (9) for each state 
in the spectrum of a ID “soft-Coulomb” potential with an 
external electric field. The parameters used are as in Fig. 
2 of Ref [16]. We used the SESM [15], with 6 = 0.4 and 
xo = 12 , lying therefore outside the simulation box (L = 10 
AU). This enables us to extract the unsealed form of the true 
resonant eigenfunctions corresponding to the perturbation of 
the lowest bound states, indicated by tpo and ipi- 


analytic derivations, makes the usage of rescaled Hamil¬ 
tonians much more powerful than the artificial addition 
of ad hoc Complex Absorbing Potentials to the device. 
This method is especially useful for coordinate mappings 
leaving the potential of the quantum device unsealed, and 
does not require any finite-difference measurement of the 
eigenvalue sensitivity, avoiding the need of tracing trajec¬ 
tories in the parameter space and to search, graphically, 
for stable points. Our findings, supported by numerical 
examples with ID model Hamiltonians, can straightfor¬ 
wardly be extended to 3D systems using e.g. separable 
forms of coordinate mappings [14], and can be applied 
to any numerical basis set. The method provides scalar 
quantities enabling the identification of physically rele¬ 
vant states and, at the same time, indicates the relia¬ 
bility of the computational setup. We believe that our 
method paves the way for accurate investigations of scat¬ 
tering properties of open quantum systems. Work is in 
progress towards the iterative extraction of Siegert states 


out of non-Hermitian Hamiltonians. 
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